#include <map>
#include <list>
#include <set>
#include <vector>
#include <cmath>
using namespace std;

#include "hla_processor.h"

class GroupsDetailProcessor : public HlaGroupsProcessor
{
public:
   GroupsDetailProcessor(){}
private:
    virtual ProcessorPriority get_priority() const  {   return 2.5f;    }
    virtual void process(const Group& group);
    virtual void generate_results(Spreadsheet& spreadsheet);

    void generate_header(Sheet& sheet);
    void generate_bottom(Sheet& sheet);

    typedef string HaplotypeName;
    typedef string SeqName;

    struct HlaData
    {
        size_t total;
        size_t unique;
    };

    typedef map<HaplotypeName, HlaData> SeqData;

    struct GroupData
    {
        string name;
        map<SeqName, SeqData> sequences;
    };

    typedef list<GroupData> Groups;

    Groups groups;
    set<HaplotypeName> haplotypes;
};

REGISTER_PROCESSOR(GroupsDetailProcessor);

void GroupsDetailProcessor::process(const Group& group)
{
    groups.push_back(GroupData());
    GroupData& gd(groups.back());

    gd.name = group.name;

    // for each haplotype in this group, analyze the peptides of each sequence:
    for (list<Haplotype>::const_iterator haplotype = group.haplotypes.begin(); haplotype != group.haplotypes.end(); ++haplotype)
    {
        // for each sequence in this haplotype, analyze its peptides
        for (list<Sequence>::const_iterator sequence = haplotype->sequences.begin(); sequence != haplotype->sequences.end(); ++sequence)
        {
            set<string> peptides;
            HlaData& hladata = gd.sequences[sequence->name][haplotype->name];

            peptides.insert(sequence->peptides.begin(), sequence->peptides.end());
            hladata.unique = peptides.size();
            hladata.total  = sequence->peptides.size();
        }

        insert_into(haplotypes, haplotype->name);
    }
}

template <class T>
static inline T sq(T x)
{
    return x*x;
}

static void calc_stats(const vector<float>& vec, float& avg, float& stdev)
{
    if (vec.size() > 0)
    {
        // first, average:
        avg = 0;
        for (size_t i=0; i < vec.size(); ++i)
            avg += vec[i];

        avg /= vec.size();

        // now, stdev:
        if (vec.size() > 1)
        {
            stdev = 0;
            for (size_t i=0; i < vec.size(); ++i)
                stdev += sq(vec[i] - avg);

            stdev = sqrt(stdev / (vec.size() - 1));
        }
        else
            stdev = 0;
    }
}

void GroupsDetailProcessor::generate_results(Spreadsheet& spreadsheet)
{
    Sheet new_sheet( spreadsheet, "Details" );
    const Style top_sep( Style::BORDER_TOP );

    const size_t n_haplotypes = haplotypes.size();
    vector<float>* const data = new vector<float>[n_haplotypes];

    generate_header(new_sheet);

    for(Groups::const_iterator group = groups.begin(); group != groups.end(); ++group)
    {
        bool first_sequence = true;
        size_t n_haplo;

        for (map<SeqName, SeqData>::const_iterator seq = group->sequences.begin(); seq != group->sequences.end(); ++seq)
        {
            Row row(new_sheet);

            if (first_sequence)
                row << top_sep << row_span(group->name, group->sequences.size()+1);
            else
                row << covered_cell;

            if (first_sequence)
                row << top_sep;
            row << seq->first << v_separator;

            n_haplo = 0;
            for (set<HaplotypeName>::const_iterator haplotype = haplotypes.begin(); haplotype != haplotypes.end(); ++haplotype)
            {
                const HlaData& hladata(find(seq->second, *haplotype));

                if (first_sequence)
                {
                    row << top_sep;
                    data[n_haplo].clear();
                }
                row << hladata.unique;

                if (first_sequence)
                    row << top_sep;
                row << hladata.total << v_separator;

                data[n_haplo++].push_back(hladata.unique);
            }

            first_sequence = false;
        }

        // output averages:
        {
            Row row(new_sheet);
            float avg, stdev;

            row << covered_cell << "Avg/Stdev" << v_separator;
            for (n_haplo = 0; n_haplo < n_haplotypes; ++n_haplo)
            {
                calc_stats(data[n_haplo], avg, stdev);
                row << avg << stdev << v_separator;
            }
        }
    }

    delete [] data;

    generate_bottom(new_sheet);
}

void GroupsDetailProcessor::generate_header(Sheet& sheet)
{
    const Style center( Style::ALIGN_CENTER );
    const Style bottom_sep( Style::BORDER_BOTTOM | Style::ALIGN_CENTER );

    {
        Row header1(sheet);

        header1 << "" << "" << v_separator;;

        for (set<HaplotypeName>::const_iterator haplotype = haplotypes.begin(); haplotype != haplotypes.end(); ++haplotype)
            header1 << center << column_span(*haplotype, 2) << v_separator;
    }

    {
        Row header2(sheet);

        header2 << "Group" << "Sequence" << v_separator;

        for (size_t haplo = 0; haplo < haplotypes.size(); ++haplo)
            header2 << center << "unique" << center << "total" << v_separator;
    }
}

void GroupsDetailProcessor::generate_bottom(Sheet& sheet)
{
    const Style top_sep( Style::BORDER_TOP );
    Row bottom(sheet);

    for (size_t haplo = 0; haplo <= haplotypes.size(); ++haplo)
        bottom << top_sep << "" << top_sep << "";
}
